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Abstract 

Background: In many biological and therapeutic contexts, it is highly desirable to target a chemical specifically to 
a particular tissue where it exerts its biological effect. In this paper, we present a simple, generic, mathematical 
model that elucidates a general method for targeting a chemical to particular tissues. The model consists of 
coupled reaction-diffusion equations to describe the evolution within the tissue of the concentrations of three 
chemical species: o (concentration of free chemical), b (binding protein) and their complex, c (chemical bound to 
binding protein). We assume that all species are free to diffuse, and that o and b undergo a reversible reaction to 
form c. In addition, the complex, c, can be broken down by a process {e.g. an enzyme in the tissue) that results in 
the release of the chemical, a, which is then free to exert its biological action. 

Results: For simplicity, we consider a one-dimensional geometry. In the special case where the rate of complex 
formation is small (compared to the diffusion timescale of the species within the tissue) the system can be solved 
analytically. This analytic solution allows us to show how the concentration of free chemical, a, in the tissue can be 
increased over the concentration of free chemical at the tissue boundary. We show that, under certain conditions, 
the maximum concentration of o can occur at the centre of the tissue, and give an upper bound on this 
maximum level. Numerical simulations are then used to determine how the behaviour of the system changes 
when the assumption of negligible complex formation rate is relaxed. 

Conclusions: We have shown, using our mathematical model, how complex degradation can potentially be 
exploited to target a chemical to a particular tissue, and how the level of the active chemical depends on factors 
such as the diffusion coefficients and degradation/production rates of each species. The biological significance of 
these results in terms of potential applications in cartilage tissue engineering and chemotherapy is discussed. In 
particular, we believe these results may be of use in determining the most promising prodrug candidates. 
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Background 

A fundamental problem in biology and medicine is how 
to target a chemical to the particular tissue where it is 
to exert its biological action. Examples of chemicals to 
be targeted include endogenous hormones, growth fac- 
tors and prescribed drugs. The advantage of targeting a 
drug to a particular tissue is that potential side-effects 
of the undesirable presence of the drug in other tissues 
are minimised. 

The results of a previous theoretical study describing 
the transport of IGF suggested a means of achieving 
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selective targeting to particular tissues. Selective degra- 
dation of the IGF-IGFBP complex was shown to be cap- 
able of regulating the free concentration of IGF (the 
biologically active form) within a tissue. In some cases 
the free IGF concentration in the tissue was raised to a 
level an order of magnitude or more greater than in the 
adjacent synovial fluid [1]. We believe that adjusting the 
rate of degradation of the complex may well provide a 
generic mechanism for tuning tissue exposure to various 
chemicals in the body. 

In order to gain further insight into these processes, in 
this paper we perform a parametric study on a generic 
system consisting of two molecules and their complex 
diffusing within a tissue. We use a combination of 
mathematical analysis and numerical simulations to 
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conduct 'mathematical experiments' to investigate the 
effect of the diffusion coefficients, degradation rates, 
binding affinities, etc. on the uptake of molecules. Con- 
sistent with the previous study of [1], the transport of 
IGF, IGFBP and their small complex in cartilage are 
considered. We also apply our model to the transport of 
a prodrug within a tumour, which provides another 
important illustration of this effect. 

The primary source of circulating IGF-1 is the liver 
and its serum concentration is relatively constant 
throughout the day in healthy adults [2,3]. The interac- 
tions between IGFs and several high-affinity IGF-binding 
proteins (IGFBPs) are important for IGF transport [4-6] 
from the liver to its site of biological action. In solution, 
binding proteins (IGF-BPs) act as 'carrier proteins', 
forming IGF-IGFBP complexes, which prolong the half- 
life of IGFs [7]. For example, in serum the majority of 
IGF-1 is bound to IGFBP3 [8]. Our recent theoretical 
study investigated the role that diffusible binding part- 
ners may play in regulating transport in a tissue [1]. 
Compared to the absence of a binding partner, it was 
shown that the existence of a complex which is broken 
down within the tissue, releasing free IGF, may lead to a 
substantial increase in the free IGF concentration within 
the tissue, even to concentrations well above those at 
the tissue boundary. Such a mechanism is biologically 
plausible as IGFBP-degrading proteases, such as matrix 
metalloproteases (MMPs), are capable of cleaving IGFBP 
into fragments that have a low binding affinity for IGFs, 
releasing functional free IGF from the complex [5,9]. 

However, there are numerous other examples where 
this mechanism for targeting a chemical to a particular 
tissue can be exploited, particularly in drug delivery e.g. 
of chemotherapeutic agents to solid tumours. To be 
effective, the drug must be able to reach the tumour in 
sufficient concentrations to kill the malignant cells, but 
it is desirable to minimise exposure of healthy cells to 
the drug, as it may cause tissue damage or lead to other 
adverse side effects. Some treatments may be targeted 
by injection directly into the tumour, but most current 
types of chemotherapy are administered intravenously, 
or orally, and so must reach their target via the circula- 
tory system [10]. This route exposes healthy and malig- 
nant cells alike to the drug. The problem is exacerbated 
by the fact that solid tumours frequently contain poorly 
vascularised regions, which means the drug must travel 
further from the blood vessel to reach its site of biologi- 
cal action. Further, advective transport with the intersti- 
tial fluid may be reduced as a consequence of leaky 
vessels within the tumour raising interstitial fluid pres- 
sure [10]. These factors, together with the fact that the 
drug is consumed within the tissue, means that those 
areas furthest from the blood vessel generally receive 
the lowest dose. 



To overcome these difficulties, a variety of approaches 
are being pursued to increase the concentration of the 
drug within the tumour, including magnetic targeting 
(which appears best suited to tumours near the surface 
of the body) [11,12], or the use of macrophages 
attracted to hypoxic regions as vehicles for drug delivery 
[13]. However, some chemotherapeutic agents take the 
form of a prodrug, which is administered in an inactive 
form and subsequently metabolised to an active form in 
the body. Experimental prodrugs (including tirapaza- 
mine (TPZ), AQ4N and PR-104) have been developed 
which are reduced to a cyctotoxic form only under 
hypoxic conditions, preventing a build-up of the active 
drug in well-oxygenated cells [10,14]. In terms of our 
model, we can view the prodrug as being a complex, 
which is broken down into an active, cyctotoxic form, 
and a by-product. 

Existing mathematical models for the transport and 
metabolism of the prodrug TPZ suggest that it does not 
penetrate efficiently into the tumour tissue [15,16], 
which appears puzzling in the light of some promising 
pre-clinical and clinical results [10] (although some 
negative results have also been reported recently [17]). 
However, these past models do not contain equations to 
describe the distribution of the active drug. Here, we 
apply our model to the transport of a generic prodrug 
within a tumour. Our analysis shows that the concentra- 
tion profiles of the prodrug and its active form may be 
dramatically different, and hence it is important to look 
at both when considering the likely effectiveness of the 
drug. This might at least partly to explain the apparent 
paradox with regard to TPZ penetration. More gener- 
ally, an improved theoretical understanding of the var- 
ious factors that affect the ability of the prodrug to 
target the tumour would allow for more rapid identifica- 
tion of the characteristics of the most promising drug 
candidates. 

This paper is organised as follows: in Methods, we 
present the model equations, which comprise a coupled 
system of three reaction-diffusion equations for the con- 
centrations of two chemicals and their complex. Model 
predictions are presented in Results and Discussion. We 
first consider the special case in which the rate of for- 
mation of the complex within the tissue is small. Under 
this assumption, we are able to solve the model analyti- 
cally, which allows us to show explicitly how the con- 
centrations of each species depend upon the various 
parameters. In this case, we are able to show that, under 
certain conditions on the parameters, the maximum 
concentration of active chemical (IGF or active drug) 
can occur in the centre of the tissue, and give an upper 
bound for this concentration. We then present numeri- 
cal simulations of the full model, showing how increas- 
ing the complex formation rate affects the solutions. 
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The paper concludes with a summary of our main 
results, and suggested directions for future research. 

Methods 

We develop a mathematical model of a simple generic 
system in which three chemicals, a, b and c, react and 
diffuse within a region of tissue, H. We assume a and b 
can undergo a reversible reaction to form a complex, c, 
the rates of association and dissociation being k+i and k. 
!. Both a and b are assumed to be removed from the 
tissue at constant rates, OL\ and a 2 , respectively. Note 
that, throughout this paper, we use the term 'removal' 
to describe all processes which result in a particular spe- 
cies being unable to participate in further interactions; 
such processes can include degradation, consumption or 
conversion to an inactive form. The complex may also 
be broken down by a reaction which results in the 
release of a, without a corresponding release of b (e.g. 
an enzymatic reaction which degrades b within the com- 
plex). This reaction occurs at a rate k_ 2 . These reactions 
are summarised diagrammatically in Figure 1. 

We denote the concentrations of the three chemicals, 
a, b and c at a position x e H by A (x, t), B (x, t) and C 
{x, t) (where t is time). Each chemical diffuses within 
the tissue, with diffusion coefficients D A , D B and D c 
respectively. Then, applying the law of mass action, we 
obtain the following system of reaction-diffusion equa- 
tions for the chemical concentrations 

9A 

— = D A V 2 A - k +l AB + (fe_! + fc_ 2 )C - a i A, (la) 
dt 

dB 

— = D B V 2 B - k+iAB - a 2 B + k- 1 C f (lb) 
dt 

dC 

— = D C V 2 C + k +i AB - (fe_i + k- 2 )C. (lc) 
dt 

Initially, we assume the chemical concentrations are 
given by 

A[x, 0) = Ai(x), B(x, 0) = Bi(x), C(x, 0) = Q(x),(2) 
whilst on the tissue boundary we set 

A = A*, B = B* f C = C* on dQ. (3) 

In terms of the cartilage example, we view a as repre- 
senting IGF, b as IGF binding protein 3 (IGFBP3), and c 
as the small, IGF-IGFBP3 binary complex. Our model 
thus represents a simplified version of models previously 
presented in [1,18]. In terms of the prodrug example, 
we take c to be the prodrug, whilst a represents the 
active form of the drug, with its removal being inter- 
preted as consumption by the tumour cells. In this 
example, b does not play an active role in the 



chemotherapy, and represents a by-product of the 
reduction of the produg. 

Nondimensionalisation 

We nondimensionalise equations (1) as follows (with 
tildes indicating dimensionless quantities) 

x = Lx, t = — , A = A*A, B = B*B, C = C*C. 
D c 

where L is typical lengthscale of the tissue, and A*, B* 
and C* are the concentrations of a, b and c, respectively, 
at the tissue boundary. The timescale chosen is the 
timescale for the complex, c, to diffuse through the 
domain. The dimensionless equations are then (drop- 
ping tildes) 

— = <5 A V 2 A - Xifi A AB + A. 2 /x A C - X 3 A (4a) 
dt 

— = 8 B V 2 B - X^bAB - X 4 B + X 5 C (4b) 
dt 

d C 

— = V 2 C + X X AB - X 2 C (4c) 
dt 

where we have introduced the dimensionless para- 
meters 

L 2 fe +i A*B* L 2 (fe_i +fe_ 2 ) 

C*D C D c 

aiL 2 a 2 L 2 L 2 \i- x C* 

X?. = , Xa. = , Xk = , 

D c D c B*D C 

„ Da D b C* C* 

&A= Dc~' &B= D- C ' ^ A = A* ' MB = F' 

The X t (i = 1, 2, 5) represent the ratios of the time- 
scale of complex diffusion to timescales of the various 
reactions: X\ is the relative rate of complex formation; 
A 2 and A 5 are the relative rates of complex breakdown 
to yield a and b, and A 3 and A 4 are the dimensionless 
rates of removal of a and b. The ratio of the diffusivities 
of a and b compared to that of c are given by S A and S B 
respectively, and (A A and [i B are the ratios of the concen- 
tration of complex to the concentration of a and b at 
the boundary. 

The solutions of the model equations are presented in 
the following section. In certain cases the model can be 
simplified sufficiently to be solved in closed form. This 
provides both validation of, and complementary insight 
to, the numerical solutions of the full model. 

Results and Discussion 

For simplicity, we consider a one-dimensional geometry, 
where the tissue is taken to occupy the region -1 < x < 
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Binding partner (b) 



OL2 



Figure 1 Diagrammatic representation of the reactions (1). 



Complex (c) 



|\> Degraded a 



Degraded £■ 



+ 



Degraded b 



1 and we assume the problem to be symmetric about x 
= 0. This geometry was assumed in previous studies of 
IGF transport in cartilage [1,18], and hence facillitates 
comparison with earlier results. In terms of the prodrug 
example, we could take it as being representative of a 
microenvironment inside a tumour where there are 
blood vessels located at x = -1 and x = 1 (see e.g. [19], 
where a similar geometry was adopted). Note that our 
assumption of symmetry about x = 0 allows us to 
restrict our attention to 0 < x < 1, with the boundary 
conditions (3) modified to 



dA _ dB _ dc _ 

dx dx dx 



at x = 0, 



(5a) 



A = B = C = 1 atx=l. 



(5b) 



We restrict our attention to the steady-state solutions 
of the governing equations, since we are interested in 
the long-term behaviour of the system. For the purposes 
of the numerical simulations, we set 



A{x, 0) = B(x,0) = C(x,0) = 1, 



(6) 



which is the simplest choice consistent with the 
boundary conditions (5). 

Analysis for small complex formation rate (k^ « 1) 

We now make the further assumption that the dimen- 
sionless rate of complex formation (Xi) is small. For the 
case of IGF, this appears plausible, since taking the 
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10' 11 M, C* 
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parameter estimates from [1] (A* 
M, L ~ 10" cm and k +1 ~ 10 5 M" 1 s" 1 ) and assuming 
to be of a similar order of magnitude to A*, and D c to 
be of the order of magnitude of D A (D A ~ 10" 7 cm 2 s" 1 , 
again from [1]) gives X\ ~ 10~ 3 . For prodrugs, we like- 
wise believe the assumption to be plausible (though it is 
not possible to give a quantitative estimate based on 
data from the literature). Typically, the prodrug must be 
converted to a radical to be active, and under well oxy- 
genated conditions, the unpaired electron on the radical 
is rapidly transferred to molecular oxygen, preventing 
build-up of the active drug. In the hypoxic environment 
of the tumour, however, this transfer cannot occur, so 
the rate at which the prodrug is regenerated is much 
lower [14]. Hence we expect the prodrug concentration 
will be much higher than those of its breakdown pro- 
ducts in the blood (so A* B*/C* will be small), whilst 
the rate of prodrug re-formation in the hypoxic environ- 
ment of the tumour will be slow (hence L 2 k +1 ID c will 
also be small). Throughout this section, we thus take X x 
« 1 in (4) (with all other parameters assumed 0(1)). 
This simplification allows us to solve the governing 
equations analytically, which provides a clear insight 
into how complex degradation affects the distribution of 
the active chemical in the tissue. The effect of relaxing 
this assumption is considered in the following section. 

For subsequent convenience, we write the solutions as 
the sum of steady state and time dependent parts, so 
that 

A = A T (x, t) + A s {x), B = B T (x, t) + B s {x), 
C = C T {x, t) + C s {x). 

Since any difference between the initial and steady 
state chemical distributions will decay exponentially fast, 
we confine our attention to the steady state solutions 
here. We expand the steady state components as a 
power series in the small parameter X\ so that 



A s = A 0 (x) + kiAi(x) + ..., 



(7) 



and similarly for B s and Cs . 

Substituting the above into the governing equations 
(4) we find the leading-order solutions obey 



a 2 A 0 



k 3 A 0 = 0, 



d 2 B 0 

8 B ^^, k 4 B 0 + A. 5 C 0 = °' 



dx 2 



dx 2 



X 2 C 0 = 0. 



(8a) 
(8b) 
(8c) 



(We note that the leading-order time-dependent solu- 
tions can also be obtained analytically by straightforward 
separation of variables and application of Sturm-Liou- 
ville theory (see e.g. [20]). The solutions of (8) are given 
by 



An 



1 + 



cosh(y 



1 + 



ti A X 2 cosh(V^), 
8 A k 2 — ^3 cosh(v^2) (9a) 



28 A 



liA\fh.x sinh(v^2^) . 



cosh(v / ^2) 



28 a cosh(V^I) 



if 8a^2 = ^3f 



n cosh 



1 + 



(5 B A 2 -A 4 )J cosh (y^ 

A.5 cosh(VXix) 
{8 B X 2 - A 4 )cosh(V^I) (9b) 



1 + 



A 5 



28 B ^k2 



tanh v k 2 



cosh +/k2~X 



cosh 



k 5 x sinh^/k^x . 



28b sfh. cosh 



if 8b^2 = 



Co 



cosh(VA2^) 



cosh(VA2) 



(9c) 



We begin by considering the case Xi = 0, for which 
the steady solutions are given exactly by (9), in order to 
investigate the effects of the other parameters. Since the 
equations for A 0 and B 0 are similar in form, we focus 
on the behaviour of the solutions for A 0 and C 0 . We 
begin by noting the following basic points. Firstly, in the 
absence of reactions {X t = 0, i = 1, 2, ... 5), A 0 = B 0 = C 0 
= 1 - i.e. the steady state solution is spatially uniform 
and the concentrations are equal to their values at the 
boundary. When there is no breakdown of the complex 
to produce a and b {i.e. X 2 = X 5 = 0), A 0 and B 0 have a 
minimum at x = 0, since they are depleted in the tissue 
due to removal. Provided X 2 > 0, C 0 has a minimum at 
x = 0, since it is only supplied at the boundary, but 
breaks down throughout the tissue. 

If the aim is to maximise the concentration of a 
throughout the tissue, intuition suggests that we need to 
reduce X 3 /S A to as small a value as possible (i.e. mini- 
mise the rate of removal relative to the rate of diffusion 
of a). Taking X 3 = 0, it is then clear from equation (9a) 
that the largest values of A 0 will occur when X 2 is large, 
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{i.e. when the rate of complex breakdown is fast com- 
pared to its rate of diffusion). Furthermore, the maxi- 
mum value of A s occurs at x = 0, in the centre of the 
tissue. In the limit X 2 — > °°, the maximum value of A 0 
— > 1+ ft A /S A , and is attained throughout most of the tis- 
sue (except for within a thin boundary layer, of width 

Q(X 2)) c l° se to x = 1, where the concentration of a 

falls to unity). Correspondingly, C 0 — > 0 in most of the 
tissue, except in the thin layer at the tissue's edge. This 
type of behaviour is illustrated in Figure 2. We note 
that, in this case, the level of a in the tissue does not 
depend upon the rate of complex breakdown (which is 
rapid), but is limited only by the supply of c to the tis- 
sue (which occurs by diffusion), and escape' of a at the 
edge of the tissue. Hence increasing the supply of c, 
either by increasing its concentration at the boundary, 
or by increasing its diffusivity results in increased levels 
of a within the tissue, as does reducing the diffusivity of 
a (which reduces the amount lost at the boundary). 
From the point of view of prodrug delivery to tumours, 
looking only at the concentration profile for c might 
lead us to think that the penetration of the prodrug is 
very poor (C 0 is only nonzero in a thin layer close to 
the boundary). However, when the concentration plot of 
the active drug (A 0 ) is considered, we would expect the 
prodrug treatment to be highly effective, as the active 
drug concentration is high throughout most of the 
tumour. 

For 0(1) values of X 3 (i.e. where the rate of consump- 
tion or removal of a is appreciable), we can observe 
three general types of behaviour. For the purposes of 
illustration, we set X 2 = 50, X 4 = X 5 = 0, \a a - ft B = S A = 
S B = 1 (recall that X\ - 0 throughout this section), and 
vary X 3 . The results are plotted in Figure 3. For rela- 
tively small values of A 3 , A s has a single maximum, in 
the centre of the tissue. (Note it is straightforward to 
show that if x - 0 is a maximum of A 0 , then A 0 has no 
other fixed points.) Conversely, for large values of A 3 , A 0 
has a single minimum at the centre of the tissue. At 
intermediate values, we see the third type of behaviour: 
A 0 has two extrema - a minimum at x = 0 and a maxi- 
mum somewhere else in the tissue. 

In terms of prodrug treatment for tumours, the aim 
would be to achieve the greatest possible level of active 
drug in the centre of the tumour (i.e. to maximise A 0 
(0)). Figure 4 shows a plot of A 0 (0) for a range of values 
of X 2 and X 3 . As might be expected, the maximum con- 
centration of a occurs when there is no consumption/ 
removal (X 3 = 0) and the rate of breakdown of the com- 
plex, c, is large. In this case we notice A 0 (0) — > 1 + fi A / 
S A = 2 - i.e. supply of c becomes the limiting factor. 

We now consider how the behaviour illustrated above 
changes when X x * 0. For small values of X\ (assuming 



all other parameters are O (1)), we can determine the 
effect of complex formation on the distributions of the 
three chemicals by considering the next-order correction 
terms A lf B x and C x in the expansion (7). The general 
solutions are 



Ai =A* cosh I / —x j + Aj cosh OjX 



+ A-j cosh(2yA^x) + As cosr^yi^c), 



(10a) 



Bi = S* cosh IJ—x) + y^Bj cosh OjX 

\ v 8b / Tt 

+£> 7 cosh(2yA^x:) + Bq + By cosh(yA^c), 



(iob) 



Ci = C* cosh yilx + Cj cosh(0j*) + C 7 cosh(2yi7x) + C 8 , (10c) 

i=i 



where the 0 t are given by: 






(Ha) 



(lib) 



(11c) 



For the sake of concision, the full details of the solu- 
tions, including coefficients Aj, Bp Cjetc. are given in the 
Appendix. 

Plotting the above correction terms for a range of 
parameter values, as illustrated in Figure 5, we observe 
that in general the effect of small but non-zero X\ is 
exactly what would be anticipated: a reduction in the 
concentration of a throughout the tissue, with a corre- 
sponding increase in the concentration of c. The effect 
is greatest in the centre of the tissue. 

Numerical simulations of the full model 

For the sake of completeness, we now turn to the gen- 
eral case where the parameter values are unrestricted, 
which requires the model equations (4) to be solved 
numerically, subject to the boundary and initial condi- 
tions (5)-(6). We use a semi-implicit finite difference 
method in Matlab. The solutions were verified by check- 
ing the results for small values of X\ against the analyti- 
cal results obtained in the previous section, and were 
found to give excellent agreement (data not shown). 
The simulations presented in this section use a timestep 
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A 2 0 0 x 

(b) Co 

Figure 2 Plots of the steady state solutions A 0 and C 0 for varying values of A 2 (see equation (9)). Other parameter values: A 3 = 0, jj A = 
100, S A = 1. Although the concentration of the complex, c, declines rapidly away from the tissue boundary, the level of a is rapidly increasing. 
Thus, for the prodrug example, apparently poor prodrug penetration into the tissue can, in fact, correspond to high levels of active drug in the 
tissue, due to efficient conversion of the prodrug to active form. 
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X 

Figure 3 Plots of the steady state solutions A 0 for various values of A 3 (see equation (9a)). Other parameter values: A 2 = 50, jj A = S A = 1. 

For small A 3 the centre of the tissue {x = 0) has the highest concentration of a; as A 3 increases, the maximum point moves towards the tissue's 
edge (x = 1), until for sufficiently large A 3 , the concentration is monotonically decreasing away from the boundary. 



dt = 0.005 and a spatial step dx = 0.002, and were run 
until the solution reached steady state. 

As we saw in the previous section, increasing Ai 
reduces the the concentration of a within the tissue, 
and correspondingly increases the level of c compared 
to the Ai = 0 case. However, when the parameter values 
are O (1), moderate increases in Ai do not make a great 
difference to the qualitative behaviour of the solutions, 
as we see from from Figure 6. It is necessary to increase 
Ai by an order of magnitude to observe such a change - 
e.g. as shown in Figure 6 for A,! = 20, it is now possible 
for A to have a second fixed point when x = 0 is a max- 
imum, in contrast to what was found in the previous 
section, for X 1 - 0. When the other parameters are not 
O (1), A and B can be large, and small values of A-i can 
make a significant difference to the concentrations e.g. 
see Figure 7 (note that in this case, our small X x analysis 
is not valid, as that requires the other parameters to be 



O (1)). But here it is important to note that the complex 
formation term introduces significant coupling between 
the solutions for A and B, so the binding partner now 
plays a role in determining the concentration of a, 
rather than simply being a passive product of complex 
breakdown. For example, repeating the simulation of 
Figure 7, but with increased X 4 (the dimensionless rate 
of removal of b) results in a decreased level of b and 
hence an increased concentration of a (compare Figures 
7 and 8). This is because there is now less b available to 
'tie up' a by forming the complex. Hence, we would 
suggest that to achieve a high concentration of a, it is 
generally desirable that the binding partner should be 
removed as rapidly as possible from the tissue. 

Conclusions 

In this study, we have considered a simple mathematical 
model for the interplay of diffusion, degradation and 



Gardiner et al. Biology Direct 201 1, 6:46 
http://www.biology-direct.eom/content/6/1/46 



Page 9 of 16 




Figure 4 Plot showing A o {0) (see equation (9a)) for a range of values of A 2 and A 3 (with \x A ld A = 1). Note that the largest values of A 0 {0) 
are achieved when A 3 = 0 (no removal/consumption of a) and A 2 is large (rapid breakdown of c), in which case A s (0) -» 1 + jJ A /S A = 2. Hence 
the highest levels of o at the centre of the tissue occur when the rate of complex breakdown (A 2 ) is rapid, and the rate of consumption of o 
(A 3 ) is minimised. 



complex formation and breakdown in the transport of 
two chemical factors and their complex. Interactions 
between these processes can produce behaviours very 
different to what is possible when complex breakdown 
does not occur. Theoretically it is clear that the key to 
targeting one of these chemicals (a) to a particular tis- 
sue (and specifically, to the centre of that tissue) is to 
ensure the complex breakdown rate is large (whilst the 
formation rate is small), and the removal rate of a is as 
small as possible. In addition the boundary concentra- 
tion and diffusivity of the complex relative to those of a 
should also be as large as possible, and the removal rate 
of the binding partner b should be large, so it does not 
persist in the tissue and 'tie up' a by forming the com- 
plex. In certain parameter regimes, it is not only possi- 
ble to raise the concentration of a in the centre of the 
tissue to a level much higher than that outside, but also 
for the concentration of a to be maximised there. This 
may be significant in terms of the treatment of hypoxic 



regions of tumours (see below). More generally, varying 
the rate of complex breakdown and the rate of removal 
can give different concentration profiles for a in the tis- 
sue. This might also conceivably be used to provide 
positional information to regulate cell behaviour e.g. 
during development. 

The results presented here complement an earlier 
modelling study which considered the transport of IGF, 
IGFBP and their complex in cartilage [1]. We have 
shown that, generally speaking, fast breakdown of the 
complex enhances the level of IGF in the tissue. How- 
ever, when the rate of breakdown becomes very fast, the 
concentration of IGF is limited by the balance of com- 
plex entering the tissue, and IGF leaving the tissue, 
embodied in the parameter ratio (Aa^a m our model. It 
is suggested that tissues may exploit the processes con- 
sidered here to 'tune' their exposure to various factors. 
As our results show, this tuning could be effected by 
adjusting the rate of complex breakdown, or the rate of 
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Figure 5 Plots showing the first-order correction terms Ai and 

Ci for various parameter values (see equation (10)). Dashed 

line: X 2 = 5, A 3 = 0.5, X 4 = 0.2, jj A = 1 . Dot-dashed line: X 2 = 2, A 3 = 

1 .5, A 4 = 0.2, ^ A = 1. Solid line: X 2 = 1 , A 3 = 0.5, A 4 = 0.2, jj A = 5. In 

all cases we take: X 4 = 0.2, X 5 = mu fi = S A = S B = 1 
v J 



removal of a (e.g. through regulation of enzyme activity 
or receptor concentration). (It appears less likely that 
the tissue would be able to influence the diffusivities or 
boundary concentrations of the chemicals, though 
changes in extracellular tissue matrix density and com- 
position may exert some effect.) Importantly, this 
mechanism can give relatively uniform levels of a 
through most of the tissue (as seen in Figure 2), which 
may be necessary for maintaining tissue homeostasis (e. 
g. for the prescence of a hormone or growth factor in a 
tissue). 

In the context of the development of prodrugs for 
treating tumours, our results suggest that not only 
might it be possible to increase the active drug concen- 
tration in the tumour to levels far above those in the 
circulation, the centre of the tumour, generally the most 
hypoxic region, can be preferentially targeted. Compared 



0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 
X 

(b) C s 

Figure 6 Plots of the steady state solutions for A and C. 

Parameter values: X^ = 0 (solid), X^ = 2 (dash), Xy = 20 (dot-dash). 
Other parameter values: X 2 = 5, X 3 = 0.5 X 4 = X 5 = jj A = jj B = S A = 
S B = 1 in all cases. Increasing the rate of complex formation {X{) 
leads to a reduction in A, and corresponding increase in C. We also 
note A can now have a maximum at x = 0 and another extremum 
in 0 <x < 1, in contrast to the X } = 0 case. 

V ) 

to some earlier models {e.g. [15,16]), which consist only 
of a reaction-diffusion equation for the concentration of 
prodrug, our model can provide new insights when 
applied to prodrug transport. Our results show that the 
distribution of the active drug may be very different 
from that of the prodrug, and that considering only the 
latter can be misleading. We have seen that low levels 
of prodrug in the centre of a tumour do not necessarily 
imply poor penetration. They may, in fact, be a very 
positive sign, associated with high levels of active drug 
in that region, due to complete conversion of the pro- 
drug to active form. This may explain the promising 
pre-clinical and clinical results of the prodrug TPZ 
referred to in [10], despite previous modelling suggest- 
ing it does not penetrate efficiently into the tumour 
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(b) 

Figure 7 Plots of the steady state solution for A and B. 

Parameter values: X } = 0 (solid), X } = 0.01 (dash). Other parameter 
values: X 2 = 1 X 3 = X A = 0.5, A 5 = 12, ^ = 10, yt B = 0.2, ^ = 0.02, S B 
= 0.1. In this case, since A and B are relatively large {B = 0(1 0 1 ) in 
both cases; data for not shown), small values of X ] can make a 
significant difference to the concentration profile. 



[15,16]. However, the rate of consumption of the active 
drug has an important role in determining its concen- 
tration profile: low rates are required for a high concen- 
tration of the active drug at the tumour centre. Higher 
rates of drug consumption drag the drug concentration 
profile down. Indeed, very high consumption rates result 
in the situation that occurs for conventional therapies, 
where the active drug concentration declines rapidly 
with distance from the blood vessel. At intermediate 
consumption rates, the drug concentration in the 
tumour is greatest at some point between the boundary 
and the centre of the tumour using this model of drug 
behaviour. Nevertheless, high consumption rates, say by 
receptor-ligand binding so dragging down drug concen- 
tration profiles, implies that the drug is actually doing 



Figure 8 Plots of the steady state solutions for A and B. 

Parameter values: X^ = 0 (solid), Xy = 0.01 (dash). Other parameter 
values: X 2 = 1 , X 3 = 0.5 X 4 = 1 0, X 5 = 12,^= 1 0, jj B = 0.2, S A = 
0.02, S B = 0.1. Compared to Fig. 7, we observe that decreasing the 
level of b (through increasing X 4 ) increases the concentration of a, 
since the rate of complex formation is reduced. 



what it is suppose to. When combined with experiments 
to determine the diffusion coefficients and reaction rates 
of the various species (along the lines of those presented 
in [15,16]), the theoretical results presented here may be 
of assistance in determining which potential new pro- 
drugs are most likely to be effective, by excluding those 
e.g. for which the diffusivity is too low, the rate of 
reduction to active form too slow, or the rate of con- 
sumption of the active drug too great. Whilst more 
complex scenarios for prodrug activation can be envi- 
saged, we believe a strength of our generic model is that 
it gives a useful insight into the fundamental principles 
involved, whilst being simple enough to be understood 
on an intuitive level. It can also provide a foundation for 
the development of more detailed models to elucidate 
more complicated situations as needed. More 
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speculatively, we believe our results suggest another 
potential avenue for drug delivery. Where a naturally 
occurring chemical and binding partner (interacting in 
the way assumed in this paper) can be found in a tissue 
of interest, a drug might be designed in such a way that 
it forms a stable complex with the binding partner. The 
drug could then piggy-back' on this native binding pro- 
tein, to be released preferentially in the tissue of interest. 

Reviewers' comments 

Reviewer's Report 1 

Dr Artem Novozhilov (nominated by Eugene V. Koonin) 
Main Comments 

In this manuscript the authors present a mathematical 
model that simulate interaction and diffusion of three 
chemicals. It is shown by means of analytical and 
numerical investigation that there are such parameter 
values that, given the initial and boundary conditions, 
it is possible to observe the maximum concentration of 
one of the chemical in the middle of the spatial 
domain, which is supposed to represent a tissue. The 
modeling results are discussed in terms of potential 
applications in cartilage tissue engineering and che- 
motherapy. In this review I would like to refrain from 
discussing possible implications of the modeling results 
with respect to, e.g., chemotherapy, and focus my cri- 
tique on the mathematical analysis of the system of 
three reaction-diffusion equations formulated in the 
manuscript. In general, I would recommend the 
authors to present their mathematical results in a 
more coherent and exact form. 
Here are my comments listed in no particular order. 

1. The authors first state that they consider the 
boundary conditions of the form A = A*, B = B*, C = 
C* (Eq. (3) on page 3). On page 4 they state that "... our 
assumption of symmetry about x = 0 allows us to 
restrict our attention to 0 < x < 1 with the boundary 

j . . , ox J<r . dA dB dc 

condition (3) modified to — = — = — = 0 at # = 0 at 

dx dx dx 

x = I" These are two different mathematical problems. 
These problems can be proved to be equivalent only in 
the case of identical initial conditions for x < 0 and x > 
0, which is generally not supposed in the paper, or in 
the case of identical asymptotical behavior of solutions 
for x < 0 and x > 0, which is not proved in the manu- 
script. Therefore the jump from the first problem to the 
second one requires some discussion. 

Response: We agree the wording of this section was not 
sufficiently precise. We have re-worded this to: and we 
assume the problem to be symmetric about x = 0/ which 
we hope makes our intention clear. 

2. The authors should elaborate on why they suppose 
that at the ends of the spatial domain all three chemi- 
cals have fixed concentrations. It can be supposed that 



one of the chemicals is targeted to the tissue; can we 
assume the same about the other two chemicals? 

Response: The assumption of constant levels of A, B 
and C on the tissue boundary, although a simplification 
is, we believe, reasonable for the following reasons. 
Firstly, in the case of IGF (and some other hormones), its 
levels in the blood remain fairly constant throughout the 
day, as mentioned in the introduction. The binding part- 
ner is also fairly constant throughout the day and so the 
complex concentrations can also be reasonably approxi- 
mated as constant. In the case of a drug treatment, we 
believe our assumption is a reasonable approximation 
for the levels in the blood when the prodrug is adminis- 
tered intravenously over a relatively long period (e.g. 
hours). It would also be reasonable for in vitro experi- 
ments where the tissue is placed in a large stirred bath 
of solution containing the treatment. Of course, in reality 
there will be examples of growth factors and drug 
administration regimes where temporal variation in the 
concentrations at the boundary is important. However, 
we wished to consider a simple, generic example here; 
such effects could be studied in future work. 

3. On page 4 the authors write (I copy the sentence 
verbatim) "We restrict our attention to the steadystate 
solutions of the governing equations, since we are inter- 
ested in the long-term behaviour of the system, and any 
different between the initial condition and the steady 
state will decay [...]" The long-term behavior of a 
dynamical system is characterized by the type of attrac- 
tor it tends to; the attractors can be steady states (both 
spatially homogeneous and non-homogeneous), cycles 
(periodic solutions) or chaotic. I don't fully understand 
why it is a priory possible to restrict the attention to 
only one type of long-term behavior. 

Response: We agree that in general several long-term 
behaviours are possible in dynamical systems, with 
steady states being just one example. However, physical 
intuition would suggest in this case the system would set- 
tle to a steady state, and in an applied paper of this nat- 
ure, we do not feel a great deal would be gained by 
establishing this by rigorous mathematical proof. We 
note our intuition is borne out by numerical simulations. 

4. The next quote from the text (page 4): "As a result, 
the initial conditions are immaterial". If there are, for 
instance, two asymptotically stable steady states then the 
initial conditions become important. 

Response: We agree with the reviewer that when there 
are two asymptotically stable steady states, initial condi- 
tions are not immaterial, and we have revised the text 
accordingly. From numerical simulations, however, we 
have not been able to find more than one asymptotically 
stable steady state in this case. 

5. The authors state (Abstract) that "In the special 
case where the rate of complex formation is small [. . . ] 
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the system can be solved analytically." Actually, the 
authors solve the system analytically only in the case 
when the rate of complex formation is zero. For this 
parameter value the system is linear, and therefore the 
steady-state solution can be written down in an explicit 
analytical form. Even in this case, however, it is neces- 
sary to prove that this steady state solution is stable to 
study the longterm behavior of the system (numerical 
experiments show that it is stable). 

Response: We agree with the reviewer notes that, in 
general when studying the long-time behaviour of dyna- 
mical systems it is important to identify which (if any) 
steady states are stable. However, considering the multi- 
disciplinary audience of this journal, we believe too 
much technical mathematical content may be off-putting 
to many readers. In the case X 2 = 0, the steady-state 
equations are linear and so their solution is unique. Of 
course, this does not prove that the system will evolve to 
this steady state, but given the simple processes involved, 
physical intuition would lead us to expect this, and it is 
borne out by the numerical simulations. We thus believe 
our approach is reasonable, although lacking mathema- 
tical rigour. We would argue similarly in the nonlinear 
(k 1 >0) case. 

6. In the case when the rate of complex formation is 
zero (Xi = 0) the authors show, using the explicit solu- 
tion, that it is possible to have the maximum of concen- 
tration at the left boundary. This conclusion is hardly 
surprising because, in words, the following system is 
considered: There is a constant inflow of the chemical 
from the right boundary, this chemical diffuses and 
degrades, and there is an impenetrable boundary at the 
left side. If the chemical does not degrade with a high 
rate it is clear that the concentration at the left bound- 
ary will increase. This conclusion does not need any 
explicit math. Therefore, it would be interesting to see 
which other insights can be obtained with the help of 
the suggested system of differential equations. 

Response: The reviewer suggests that a verbal argu- 
ment suffices to understand the fact that, when X x - 0, 
the active drug concentration can have a maximum at 
the centre of the tissue. We agree, but would point out 
that not all the results of the model, even in this simpler 
case, may not be so clear on an intuitive level. For exam- 
ple, the fact that the maximum drug concentration may 
occur at a point intermediate between the centre and the 
edge of the tissue. If such an observation was made 
experimentally, we do not believe it would necessarily be 
intuitively obvious that it could be produced by the sim- 
ple mechanism outlined here. We also see there is a 
maximum theoretically attainable active drug concentra- 
tion (A* + D c C*/D A in dimensional terms), which again 
might not be intuitively clear, and obviously has practi- 
cal implications. 



7. The authors argue that the found analytical solution 
is very close, in case of small X lf to the steadystate solu- 
tion observed in the numerical experiments. The agree- 
ment, according to the text, is excellent, but they do not 
show the data. Let me consider the following values of 
the parameters: Xi - 0.1, X 2 = 1, X 3 = 0.5, A 4 = 0.1, X 5 = 
12, [A A = 10, ft B = 0, S A = 0.01, S B = 0.1. Using the exact 
formula (9a) from the text and the standard numerical 
PDE solver, I obtain the following figure (Figure 9). 

Here the blue line shows the exact solution for A 0 (x) 
of the corresponding linear system, the red line shows a 
numerical solution for A (x) of the full nonlinear system 
for large t. This figure gives a clear example that it is 
dangerous to make any conclusions about the long-term 
behavior of the nonlinear system using the exact solu- 
tions of the approximate linear system. From a mathe- 
matical standpoint it is only possible to state that the 
solutions of two systems are e-close on the times of 
order 0(1 /e) but nothing can be proved for t — > °°. In 
general, the authors should either not consider the exact 
solution of the approximate linear system or specify the 
parameter domain in which the linear approximation 
works. 

Response: We thank the reviewer for drawing our 
attention to this point. As is stated in the section on the 
analysis of the small X x case, we need to assume that the 
other parameters are 0(1) for analysis to be valid, and 
we have now reiterated this at the beginning of the 
appendix. (In passing we remark that in some cases, the 
approximate, small-Xi solution gives good agreement 
even when this assumption is violated; however, this can- 
not be expected in general.) We note that in the case pre- 
sented by the reviewer, the ratio [t A /d A (which appears in 
the equation for A) is O(10 3 ); hence the grouping X x fl A / 
S A is not small, and we would not expect our analysis to 
be valid. However, we have now included an example of 
this type of behaviour (the new Figure 7) in the numeri- 
cal results section, as it demonstrates the important 
point that there are situations where a small value of X x 
can nevertheless make a significant difference to the solu- 
tion. For the cases where the other parameters are 0(1), 
we have compared our analysis with the numerics for 
several cases. The agreement with the leading-order solu- 
tion is generally reasonable, and improves with the 
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Figure 9 Numerical simulation provided by Reviewer 1. 
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addition of the first order correction term. However, since 
this does not really provide any new information to the 
reader, the relevant plots are omitted. 
Minor Comments 

1. Page 2, second column. . . the drug must the. . . " 
should be ". . . the drug must be." 

2. Page 3, second column. Why use indexes for the 
initial conditions? 

3. Page 5, second column. . . exactly by (8). . . 
should be ". . . exactly by (9)." 

4. Page 6, second column. A missing number in . . 
in the expansion ()." 

5. The first reference I checked has a typo: the paper 
by Zhang et al. was published in 2010, not in 2009. 
Please check all the other references. 

Response: We thank the reviewer for also pointing out 
a number of smaller errors, which we have corrected 
according to their suggestions, except for retaining the 
subscript i for initial conditions, where we think the 
intended meaning is clear. 

Reviewer's Report 2 

Dr William Hlavacek 
Main Comments 

The authors present steady-state results, including ana- 
lytical results, for a model consisting of three partial dif- 
ferential equations for a system in which there is a 
reversible bimolecular association reaction (A+B = C) 
and one-dimensional diffusion of the chemical species 
A, B and C. The analytical results are approximate and 
obtained through a perturbative approach. The analyti- 
cal expressions given by the authors are fairly compli- 
cated, and it is unclear that these expressions provide 
more insight than purely numerical results. The model 
is simple, but the authors explain that models even sim- 
pler in some respects have been considered in earlier 
related work and that their study is relevant for under- 
standing penetration of the activated form of a prodrug 
into a solid tumor and for understanding growth factor 
action in tissues. The report is well written. The insights 
gained from analysis of the model seem somewhat 
obvious. For example, the authors conclude from their 
analysis that rapid clearance of the activated form of a 
prodrug will decrease its effectiveness. I'm not sure any- 
one would need a model to understand this point. Also, 
it is unclear that any of the insights gained from the 
analysis are actionable. For example, the authors suggest 
a piggy-backing strategy for increasing drug penetration 
but no concrete example is provided. Can the authors 
provide more details about potential proteins that could 
serve as carriers of drugs into solid tumors? In any case, 
the results are purely theoretical. There is no demon- 
stration that the results are useful in practice. There is 
only limited discussion of parameter estimates. 



Response: We thank the reviewer for their report. They 
note that the fact that e.g. increasing the rate of active 
drug clearance tends to reduce its concentration in the 
tissue is intuitively clear. We agree, but would point out 
that there are phenomena which can be observed in the 
model which may not be so clear on an intuitive level. 
For example, the fact that the maximum drug concentra- 
tion may occur at a point intermediate between the cen- 
tre and the edge of the tissue. If such an observation was 
made experimentally, we do not believe it would necessa- 
rily be intuitively obvious that it could be produced by 
the simple mechanism outlined here. Of course, many 
other processes may produce such concentration profiles, 
and each case must be considered individually. However, 
we believe in general it would be sensible to consider the 
possibility of a scenario of the type we have suggested. 
We also agree that the analytical expressions (which are 
given for the sake of completeness) are somewhat compli- 
cated, and in are best illustrated by plots as for the 
numerical results. However, we believe they can give use- 
ful information about the dependence of the solution on 
the various parameters. For example, we see there is a 
maximum theoretically attainable active drug concentra- 
tion (A* + D c C*/D A in dimensional terms). This is 
potentially useful practically, since if this level is below 
that at which the drug is effective, alternative strategies 
for drug delivery will need to be investigated. Similarly, 
our model suggests a framework by which simple in vitro 
experiments might be used to characterise which are the 
most promising drug candidates for development by 
measuring parameters such as their diffusion coefficients 
and reaction rates: e.g. those with the smallest value of 
Xi or A 3 or largest X 2 appear more likely to give rise to 
effective treatments. We note the recent study announ- 
cing negative results in a trial of TPZ in treating 
advanced squamous cell carcinoma of the head and 
neck, and have included a reference to this paper. How- 
ever, the argument we wished to make in our paper was 
a general one: that considering merely the concentration 
of the prodrug (C in our model) could be misleading, 
since the concentration profile of the active drug (A) in 
the tissue could be very different. This applies to any 
drug which is delivered in the manner considered here. 
Owing to lack of experimental data, we cannot make 
any definitive statements about the effectiveness of TPZ, 
and the reference to it was merely offered as a possible 
example. We have re-written this paragraph of the 
'Background' section to try to make our intention clearer. 

Reviewer's Report 3 

Dr James Faeder 
Main Comments 

This paper presents a simple mathematical model to 
describe the spatial distribution of drug in tissue where 
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the drug may interact with another molecule (e.g. pro- 
tein) to form a complex that undergoes degradation at a 
different rate from the drug alone, producing back again 
the active drug and eliminating the binding partner. The 
appropriate partial differential equations to describe this 
system in the continuum limit are derived and solved 
exactly or approximately in limiting cases, and these 
solutions are supplemented by the use of numerical 
simulation on the full model. The main findings are that 
under appropriate conditions, e.g., when the degradation 
of active drug is relatively slow, the breakdown of the 
complex is fast, and the rate of complex formation is 
slow relative to the diffusion constant of the drug, the 
concentration profile of the drug in tissue reaches a 
maximum at the point farthest away from the bound- 
aries (under the assumption of fixed concentrations on 
the boundaries and a symmetric, ID tissue). Such a pro- 
file suggests a potential novel mechanism for the activity 
of a drug that must be targeted selectively to a tissue by 
utilizing enhanced tissue-specific degradation of a com- 
plex between the drug and another molecule or protein. 
Alternatively, some other tissue-specific factor, such as 
hypoxia in the case of a tumor, may be used for rate- 
enhancement of the activating process. I found this 
paper to be very clearly written both at a general and 
technical level. The findings are also clearly significant, 
both for their potential therapeutic implications and for 
modelers interested in applying this model or extensions 
of it to model other systems. I have major criticisms of 
the paper. - At the bottom of 2nd column on p. 3., I 
found the sentence "In this example, b does not play an 
active role. . . " to be confusing. I think it would help to 
explain the application of the model to the prodrug 
example a bit more clearly. Specifically, how is b mod- 
eled in that case? 

-There is a typo in Eq. 8c, which is missing '=0' before 
the period. 

Response: We thank the reviewer for their positive 
report We have addressed the minor matters raised in 
the review. 

A Details of the solutions for A lr and C-, given 
in equation (10) 

For the sake of completeness, we give here the details of 
the solutions for the first-order correction terms (note 
that we assume that the other parameters, apart from Ai 
are 0(1)). Substituting the expansions (7) into the gov- 
erning equations (4), at 0(Ai), we find 



d 2 Ai ULa ^2t^A ^3 

—± - ^A 0 B 0 + -j^Ci - /Aj = 0, 

OX z OA Oa Oa 

d 2 B\ lir Xa _ 

-A 0 B 0 - -ifi! + /Ci = 0, 



dx 2 8g 

dx 2 



+ A 0 B 0 — X 2 Ci. 



The above are to be solved subject to the boundary 
conditions 



0 at x = 0, 



3Ai _ dBi _ dCi 

dx dx dx 

A x =Bi = Ci = 0 atx= 1. 



The general form of the solution to the above is given 
by (10), whilst the details of the coefficients are provided 
below. The Aj and A* are 



A r - 



A 7 - 



, A 



for i = 1, . . . , 6, 

-fi A x 2 c* 



cosh t 



Aj cosh 6j + A-] cosh l^/Xj + As coshf-v/^e) 



C V =1 



where the 6j are given by equation (11), and the coeffi- 
cients Cj and C* are given below. 
The coefficients in the expression for B x are given by 



B 2 = - 



B 4 = - 



B 7 - 





8 B (el - 


x 4 \ 
s B ) 


s B (el - 


xA 
s B ) 


8 B (el - 


x 4 \ 
s B ) 


8 B (el - 


x 4 \ 
8b ) 


8b (el - 


x 4 \ 
8b ) 


8b (el - 


X 4 \ 
Ss ) 


8 B (4A 2 - 
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8b ) 







( HbPiPa 



-A. s Ci 
-k 5 C 2 

- A.5C4 
-k 5 C 5 

- ^5^6 

-k 5 C : 



R 1 { -HbPiP* r \ 



Cs), B 9 



k 5 C* 



8 » 



cosh 7l V 
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where we have introduced 



Pi = 

Pi 



cosh J ^ 



1 + 



1 



^2 Ma 

8a^2 — ^3 
X 2 [I A 



COSh V^2 V <$A^2 — ^3 7 ' 

1 A A 5 




04 = " 



COSh V^2 V $B^2 — ^4 , 

Finally, the coefficients in the solution for Ci are 



C 5 



2(0? -* 2 )' 
Pi Pa 



2(ei -x 2 y 



c 7 



PiPa 
' 6X 2 



C 2 - 
Ci- 
Ce- 
Cs 



2{ei-X 2 )' 

Pi Pa 
2(61 -X 2 )' 



2(6i-k 2 ) 
Pi Pa 



2X ? 



cosh 



^2 Cj cosh Oj + C 7 cosh 2^1^ + C 8 



We note, as is clear from the above, that the solution 
(10) is only valid provided Of i ^etc; however, the solu- 
tions in the special cases where this doesn not apply can 
be found from (10) by straightforward application of 
L'HopitaFs rule. 
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